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Abstract 

We simulate the dc magnetic response of the diluted dipolar-coupled Ising magnet 
L1H00.045Y0.955F4 in a transverse field, using exact diagonalization of a two-spin Hamiltonian aver- 
aged over nearest-neighbour configurations. The pairwise model, incorporating hyperfine interac- 
tions, accounts for the observed drop-off in the longitudinal (c-axis) susceptibility with increasing 
transverse field; with the inclusion of a small tilt in the transverse field, it also accounts for the 
behavior of the off-diagonal magnetic susceptibility. The hyperfine interactions do not appear to 
lead to qualitative changes in the pair susecptibilities, although they do renormalize the crossover 
fields between different regimes. Comparison with experiment indicates that antiferromagnetic 
correlations are more important than anticipated based on simple pair statistics and our first- 
principles calculations of the pair response. This means that larger clusters will be needed for a 
full description of the reduction in the diagonal response at small transverse fields. 
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I. INTRODUCTION 



The dipolar rare-earth magnetic salt LiHoF4 orders at 1.53 K-^ to form a Ising-like 
ferromagnet with long, needle-shaped domains oriented along the Ising c-axis 3 -. This system, 
and the dilution series LiHo x Y 1 _ :r F4 with the magnetic Ho 3+ ions replaced by non-magnetic 
Y 3+ , have been studied for more than three decades^ - — . For moderate dilution (x > 30%) the 
system continues to behave as an Ising ferromagnetic^; however for smaller x it appears 
to form a spin glass at low temperatures^. At x = 4.5% there is evidence for a novel 
antiglass^- in which the scaled distribution of relaxation times loses its low-frequency tail as 
the sample cools. In this phase the material exhibits macroscopically long-lived magnetic 
excitations^ and a novel combination of strong features in the specific heat with a featureless 
magnetic susceptibility which can only explained by positing long-range spin entanglement^. 
Some recent experiments have reported contrasting results — notably a featureless specific 
heat from x = 1.8% to x — 8%r 3 -, suggesting that the conventional spin glass may persist to 
lower concentrations. 

The dynamics in these dilute phases are particularly interesting and could well be the 
key to understanding the seemingly contradictory experiments. As well as the long-lived 
magnetic oscillations revealed by hole-burning experiments at x — 4.5%^, cotunnelling of 
the electronic and nuclear moments on pairs of neighboring Ho 3+ ions has been observed at 
x = 0.1 %i£ through its effect on the low- frequency zero- field susceptibility. It is appropriate 
to revisit the low-frequency susceptibility for several reasons. First, LiHo x .Yi_ x F 4 is expected 
to be a model for a wide class of transverse-field dipolar systems. Second, the observation 
of long decoherence times and signatures of long-range entanglement suggest the possibility 
of exploiting the Ho 3+ ions as magnetic qubits. Finally, one would like to understand the 
precise role of the competition between the collective dipolar interaction, the nuclear spin 
bath and other decoherence pathways in determining the dynamics of the system 19 . Here we 
combine an experimental study of the magnetic response of the dilute system as we tilt the 
moment away from the Ising axis under large transverse fields with a theoretical analysis in 
which we average over all possible pairs. Our purpose is to establish — quantitatively — the 
extend to which collective (i.e. beyond-pair) effects are important for the behavior of the 
x = 4.5% compound by doing the mos precise possible calculations of the pair susceptibility 
contribution at equilibrium. The outcome is that even for this relatively high level of dilution, 
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the collective effects are important at low transverse fields. 

We presented the experimental results and a short summary of the theoretical argument 
This paper gives full details of the model and is structured as follows. Section |TT] 



in Ref . 



18 



summarizes the experimental techniques employed and captures briefly the relevant results; 
Section [TTT] describes the techniques employed in our calculations; Section [TV] sets out the 
computational results, comparing the susceptibilities with and without hyperfine interactions 
and comparing them to the measured values; and Section |V] presents our conclusions. 



II. SUSCEPTIBILITY MEASUREMENTS 

A single (5 x 5 x 10) mm 3 crystal of LiHoo.o45Yo.955F4 was characterized using ac magnetic 
susceptibility in a helium dilution refrigerator. The magnetic response along the Ising axis 
and in the transverse plane were experimentally measured using a specially devised multi- 
axis ac susceptometer, as shown in Fig. [TJ The sample was probed using a 101 Hz 2 /iT ac 
magnetic field parallel to the Ising axis. A pair of nested inductive pickup coils allowed for 
simultaneous determination of the magnetic response parallel to and transverse to the Ising 
axis of the crystal. The crystal was thermally linked to the cold finger of the refrigerator 
via sapphire rods and heavy copper wires. A multi-axis set of 100 mT Helmholtz coils and 
an 8 T solenoid provided dc magnetic fields Hd c parallel to and almost transverse to the 
Ising axis respectively; however, because of the difficulty in precisely aligning the crystal, 
we cannot exclude the possibility of a small (~ 0.6°) misalignment of the solenoid from the 
transverse axis. The effects of such a misalignment on the predicted properties are discussed 
in CVDl below. 

The measurement probes the diagonal and off-diagonal components respectively of the 
linear susceptibility tensor, but evaluated at the non-zero reference field H dc : 

dM.. 



Xz 



x., 



dH z 
dM r 



dH z 



; (i) 

H— H(j c 

(2) 

H— H dc 

(3) 



Fig [2] shows our results for the real part of the longitudinal and transverse susceptibilities 
Xzz an d Xxz as functions of Hd c 18 . These experimental results will be compared in section 
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FIG. 1: (Color online) Schematic of the ac vector susceptometer used in the experiments. The 
sample sits inside nested pickup coils A and B, sensitive to magnetic response in the transverse 
and Ising directions, respectively. An ac magnetic field along the Ising axis is supplied by solenoid 
C; the sample is thermally sunk to the cryostat cold finger via sapphire rods D and copper wires 
E. A superconducting 3- axis Helmholtz coil F and an 8T solenoid magnet G supply dc magnetic 
fields. G is almost, but not perfectly, aligned transverse to the c-axis of the sample. 

HVDI to the predictions derived from the spin-pair model developed in the following sections. 
The off-diagonal linear susceptibility vanishes in the limit where Hd c is exactly perpendicular 
to the Ising axis; as we shall see, a small component along z enables \ xz to capture some of 
the non-linear dependence of M on H and hence to give information about clustering and 
correlation effects, as expected from previous work—. 

The imaginary part of the magnetic response was also measured. Since the frequencies 
involved are small compared with all the energy scales of the microscopic Hamiltonian, a 
theoretical treatment of the dissipation depends on an understanding of the low-frequency 
relaxation dynamics of the Ho 3+ ions and is not considered in the present paper. 
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FIG. 2: (Color online) Measured longitudinal (top) and transverse (bottom) real susceptibility at 
70, 110, and 150 mK (blue). (Adapted from Ref.iS) 

III. HO s+ PAIR MODEL 

To construct a model for the susceptibility of Ho 3+ pairs, we start with the complete 
microscopic Hamiltonian. The low-lying states of this Hamiltonian are then used to construct 
an effective 2-state H, which can be readily diagonalized for two interacting ions. If the 
hyperfine interactions from the microscopic single-ion Hamiltonian are added to this 2- 
state picture, the resulting H has 16 states, and the pair Hamiltonians are still numerically 
tractable. Finally, a weighting scheme is implemented that incorporates contributions for 
pairs beyond immediate nearest-neighbors. 

A. Microscopic Hamiltonian 

The electronic Hamiltonian of a single Ho 3+ ion in a magnetic field is 

Hi = H ci — m ■ B 

(4) 

= H ci - /iBfi'ij • B , 
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FIG. 3: Single-ion energy levels as a function of longitudinal magnetic field, (a) Lowest eight 
electronic crystal-field levels of the 5 Is ground term as a function of field B z parallel to the Ising 
axis, (b) Splitting of the two lowest electronic levels by the hyperfine interaction. 



where = | is the Lande g factor. H c f is the crystal field Hamiltonian, which splits the 



i 

5; 



17- fold degenerate J 8 ground term state of Ho, and is given by 

H ci = £ B i°i + E B t^)OKc) + Bf(s)Ot(s) , (5) 



l=2A,6 1=4,6 



where Of 1 are Stevens' operators . We follow Ref. 
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in taking the following values for the 



crystal-field parameters: 5° = —0.06 meV, = 3.5 x 10 4 meV, B\ = 3.6 x 10 3 meV, 
Bl = 4 x 10~ 7 meV, 5 4 (c) = 7.0 x 10~ 5 meV and £ 4 (s) = 9.8 x 10~ 6 meV. The resulting 
electronic energy levels are shown in Figs. |3£i and as a function of fields parallel and 
transverse to the Ising axis. 

The isotropic hyperfine coupling to the local / = | Ho 3+ nuclear spin can be included 
explicitly by defining 

H M = H ct <g> I N + A3 ■ I + fi B g L J • B + // N I • B, (6) 

with J a = j a (g) I N and A/k B = 0.039 K or A = 3Afj,eV. Figs. [3b and Hb show the effect 
of the hyperfine splitting on the lowest two crystal-field states (but computed using the 
entire single-ion Hamiltonian (jHJ)). As emphasized by Ronnow et. al. 19 and Schechter and 
Stampi 1 -^, although A is small compared with the characteristic intra-ion electronic energy 



scales, it is comparable to the inter- ion dipolar coupling (see §IV CI) . Its effect is to suppress 
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FIG. 4: Single- ion energy levels as a function of transverse magnetic field, (a) Lowest three 
electronic crystal-field levels in the presence of a field B x transverse to the Ising axis, (b) Splitting 
of the two lowest electronic levels by the hyperfine interaction. 

the mixings between the two terms of the lowest electronic doublet at low temperatures, 
because the lowest electro-nuclear spin state in each branch has the nuclear and electronic 
moments anti-aligned and the nuclear moments cannot be reversed at low orders by any of 
the terms in equation (J6j). 

The state-space required to correctly describe the 5 Is ground term of Ho 3+ in the presence 
of hyperfine splitting is then (2 x 8 + 1) x (2 x | + 1) = 136. The full Hilbert space on 
an ion pair therefore has dimensionality 136 2 = 18496, which is inconveniently large for 
the repeated exact diagonalizations required to treat a range of pair geometries and fields. 
We therefore proceed by truncating the model to a smaller state space while preserving the 
essential behavior. 

B. The electronic two-state system 

Following Chakraborty et al. 2 , we note the large (9.5 K) gap between the ground state 
doublet and the first excited crystal-field level (FigHk). We therefore construct a Hamilto- 
nian describing the low-energy behavior of the ion on a two-dimensional electronic Hilbert 
space, covering only these states. This is a parameterized model in which the inter-level 
repulsion shown in Fig. is included explicitly as described below. 
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For a given value of transverse field B x the following two-state Hamiltonian is defined: 



Here I2 is the identity operator in two dimensions and a x is a spin-half Pauli operator. 
Eq(B x ) is the mid-point of the lowest two energy levels and A(B X ) their splitting in that 
transverse field. The effective angular momentum operators are chosen to reproduce the 
correct physical angular momentum matrix elements for the two states; their decomposition 
into Pauli operators is discussed in Ref. (2]. Finally, the field B has been replaced with 
B' = B - B x i. 

Note that at first sight one might expect that it would also be possible to construct a 
three-state model, including the two-fold degenerate ground state as well as the first excited 
state, which are relatively well separated from the rest of the spectrum (see Figure [3]). 
However it turns out that level repulsion from the rest of the spectrum becomes significant 
at modest external fields-, and for this reason it is preferable to parameterize a two-state 
effective Hamiltonian operator for every value of transverse field in order to incorporate all 
these effects. 

In the presence of the I = ~ hyperfine interaction, the two-state model becomes 



with J e ff = j e ff ® In. This has a dimensionality of 16, and thus the Hamiltonian of a pair 
of spins will have a numerically tractable dimensionality of 256. In this paper we therefore 
retain the full nuclear Hilbert space when considering the hyperfine interaction, rather than 
restricting the model further to the lowest electro-nuclear doublet as in Ref. I22J. 

C. Intra-ion coupling 

We neglect the small exchange interactions between the Ho 3+ ions, so in our model pairs 
are coupled only by the magnetic dipole interaction. Angular momentum operators are 
constructed for each spin in a direct product Hilbert space. The dipole coupling between 
spins at Rf and R2 is then 





H$ = E (B x )I l6 + -A(B x )a x <g> I 8 + a^iAs ■ B' + /x N I • B + AJ oS I , (8) 





S 



where R42 = R2 — R-i and J a is component a of the total angular momentum of ion i. The 
total Hamiltonian of the pair is 

tfpair = H 1 + H 2 + H 12 . (10) 

Note that for a given pair, if 12 gives rise to an effective field at site 1 



pi J) M/^BffL) ( X 3 ^12^12 \ ,<->> ...... 



which in general contains a transverse component. Strictly, therefore, the field-dependent 
parameters in equation (0) should be computed incorporating this component. However, in 
practice this dependence is negligible for the applied fields of interest because the character- 
istic scale of is at most ^o^b9l\J^\/o- 3 = 29 mT, while the experimental variation of x 
is on the scale of fields that can mix the Ising doublet, of order IT (see Figures [2] and H]) . 



D. Computing the susceptibility 

The isothermal susceptibility is defined as 

_ 1 ( d(m Q ) \ 

Xaf3 = v{~dWi) T > (12) 

where m is the total magnetic moment and V is the sample volume. We apply this by 
computing the field-dependent eigenstates of the pair Hamiltonian ( flOl) and computing 

X«/3 = - kB Tzv ^ ex P(~^/ k B T ) (i|Am a |i> (i\Arhp\i) 



Ei — Ej 



(13) 



ZV 

i 

= XPauli + XVan Vleck , 

where the primed sum goes over all states % and j such that Ei 7^ Ej and Am a = m a — (m a ). 
Matrix elements between degenerate states have been made to vanish by a choice of basis 
such that rhp is diagonal in each degenerate subspace. Numerically we assume states % 
and j are degenerate if Ej — Ei < e, a small value chosen such that the susceptibility is 
not sensitive to variations in e. Note that in applying equation (1121) we assume that the 
Ho 3+ ions remain in thermal equilibrium over the timescales of the experiment, i.e. that all 
thermalizing relaxation processes operate on a timescale fast compared to the measurement. 



E. A pair-ensemble weighting scheme 



We wish not only to examine the behavior of specific pairs of spins, but also to cal- 
culate the average response for a distribution of spin pairs corresponding to the physical 
LiHo .o45Y .955F 4 crystal. We proceed by assuming that at this dilute concentration the 
behavior of each spin is affected only by the closest spin and compute a weighted suscepti- 
bility. This is determined by computing the susceptibility of an exhaustive sample of pairs 
of spins up to some cutoff distance r c and weighting each term by the probability that in a 
randomly populated set of sites in a lattice with mean fractional occupancy x, the chosen 
spin s 2 would be the nearest occupied site to the reference spin si. If all the sites were at 
different distances, this would be given by the probability that no sites nearer to Sj than s 2 
are occupied, while the site s 2 itself is occupied. The weighting for a site Sj would then be 

Wj = x{l - x) N * , (14) 

where Nj is the number of sites closer to s\ than sj. However in practice the sites s 2 occur 
in 'shells' with equal distance from si; if there are rij sites in shell j, we ascribe a weighting 
to each site which is a fraction 1 / rij of the probability that there is at least one neighboring 
spin anywhere in the shell: 

(1 (15) 

The cutoff distance r c is always chosen such that the probability of the nearest occupied 
site s 2 being more than r c from si does not significantly exceed 10 -3 ; the required r c therefore 
increases as x falls. For the calculations presented here we included 22 shells of neighbors 
containing 146 ions, corresponding to r c = 2.58 a = 13.4 A. At the experimental spin 
concentration (x = 0.045) the probability that the pair separation exceeds r c is then 1.20 x 
lO- 3 . 

IV. RESULTS 

A. Contributions of individual pairs 

The magnetic response of a pair of Ho spins depends strongly on their separation and 
orientation. Fig [5] shows the Ising-axis and transverse response of all pairs that make a 
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(1 - xp 



significant contribution to the cluster ensemble. Although these plots are of illustrative 
value in demonstrating the wide range of behaviors arising from spin pairs, it is more useful 
to examine how these different responses contribute to the ensemble average. Fig shows 
these averages by plotting the susceptibilities of each pair using the weighting Wi as a color 
map. Susceptibility bands appear in this weighted map due to particular closely neighboring 
spin pairs. It can also be seen that for every pair of spins with a transverse susceptibility 
Xxz{B x ) = f(B x ), there exists a pair with Xxz{B x ) = —f(B x ). It thus follows that an 
ensemble average as defined in Sec. IIIIEI will give a zero value of Xxz for all values of field 
B x . As discussed below, the measured response is well described by a small (0.6°) tilt 
of B x , producing a polarizing field along the Ising axis. A comparison of the weighted 
susceptibilities with and without the incorporation of hyperfine effects suggests that the 
primary effect of the hyperfine term is to renormalize the transverse field; this behavior is 
discussed in more detail in Section IIVCI below. 

B. Pair orientation and response 

Depending on relative orientation, the dipole coupling can be either ferromagnetic or 
antiferromagnetic. A ferromagnetic pair has a susceptibility Xzz which diverges in the limit 
of low temperatures and zero transverse field, whereas an antiferromagnetically coupled pair 
has vanishing susceptibility in the same limit. As can be seen from Fig. [2} antiferromagnetic 
behavior dictates the measured response of the sample of LirIo .o45Yo.955F 4 , and as shown 
in Fig [5] certain pairs show a qualitatively similar magnetic response. As we shall see below, 
however, their contribution to the ensemble average used in this paper is not sufficient to 
make the overall average susceptibility agree with the measured one. 

The relation of this behavior to the crystal geometry can be understood from Fig. 
showing the zero-field susceptibility at T = 70 mK of a pair of Ho 3+ ions separated by a 
distance r in the a-b plane and z on the c-axis. The crossover between the ferromagnetic and 
antiferromagnetic couplings occurs along the line z/r = l/v2; the strongly antiferromag- 
netic pairs are located in-plane at (1, 0, 0) and (2, 0, 0) and the most strongly ferromagnetic 
pair is the nearest-neighbor pair at (1,0, 4). Note that the on-axis pair (0,0,1) is more 
weakly ferromagnetic at this temperature, owing to the larger spatial separation. 
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FIG. 5: Computed susceptibilities for all pairs of spins at T=70 mK, with hyperfine interactions 
included, (a) Diagonal response Xzz- (b) Off-diagonal response Xxz- 

C. The effect of the hyperfine interaction 

We now examine the role that hyperfine interactions play in determining the behavior of 
the system. It is important to understand whether these effects produce a qualitative change 
in the behavior, as expansion of this model to n = 3 and larger clusters of spins becomes 
numerically impractical if the hyperfine splittings are essential. Fig. [8] shows susceptibilities 
for high-weight spin pairs both with and without hyperfine effects. (Note that pairs such 
as (|, 0, |) and (0, ~, |), which are equivalent at zero field, become inequivalent for non-zero 
fields, except when the field lies along symmetry directions such as (1,1,0).) We see that 
the primary role of the hyperfine interactions is to renormalize the applied transverse field, 
rather than to introduce fundamentally different behavior. This in turn suggests that useful 
insights may be derived from considering larger spin clusters in the absence of the hyperfine 
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FIG. 6: (Color online) Contribution of the various pairs to the ensemble-averaged functions Xxz(B x ) 
(top) and Xzz(B x ) (bottom). Left and right, respectively, show the effects of omitting and including 
the hyperfine term in the Hamiltonian. The temperature was T = 70 mk and the field was applied 
along (1,0,0). 

splittings. It should be noted, however, that the strongly ferromagnetic (n,0, 4) pair does 
not show this renormalization when it is oriented so that the projection of the separation 
vector into the a6-plane lies along the transverse field direction. 

D. The ensemble-averaged susceptibilities 

Fig. [9]shows the experimental and ensemble-averaged longitudinal susceptibility Xzz- The 
left panel shows computed and experimental results at temperatures of 70, 110 and 150 mK. 
Computed results include the effect of the hyperfine response, but omit in this panel the 
effect of tilting the field B x . The model captures the overall temperature dependence of the 
data, but it cannot account for the low-field suppression of the susceptibility because the 
average is dominated by the contributions of ferromagnetic and effectively uncoupled pairs. 

The right panel of Fig. [9] shows the effects of varying the parameters of the model at a 
constant T=70 mK. The dashed curves show the result of removing the hyperfine terms; 
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FIG. 7: (Color online) The effect of geometry on Xzz, computed with zero transverse field and 
temperature T = 70 mK. Response is plotted for a pair of spins with axial separation z and in- 
plane seperation r (units of lattice parameter a) with the marked points showing the locations of 
various nearest neighbors. The susceptibility is shown in units of emu/mol Ho. 

for most of the field range, the renormalization seen in the individual pair susceptibilities is 
visible. At low field, the strongly ferromagnetic (^,0, |) pairing dominates, and no renor- 
malization is seen. The dotted curve shows the result of keeping the hyperfine effects and 
adding a 0.6° tilt to the applied field, with the attendant slight polarization along the Ising 
axis. We can see that this improves the match between the high-field behavior of the model 
and the experiment. 

Fig. [10] displays similar information for \ xz . Note that owing to the symmetry observed 
in Fig. Mjo), the ensemble average of \ xz vanishes in the absence of a polarizing field. Thus, 
the only appropriate comparison is between the tilted-field computation and the measured 
value, as shown in the left pane of Fig. [10] for both single-ion and ensemble-pair-average 
computations. It is clear that the tilt is responsible for the measured effect, with the pair 
average providing a better match to the measured susceptibility than a single-ion calculation. 
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FIG. 8: (Color online) The magnetic response at T = 70 mK of certain important spin pairs, using 
a Hamiltonian which incorporates hyperfine effects (solid) and which omits these effects (dotted). 
The primary effect of adding the hyperfine splitting is to impose an effective renormalization of 
the transverse field scale. The transverse field is applied along the (1,0,0) direction. 

The effect of the hyperfine response is the same renormalization of the field seen in the 
longitudinal response. The right pane of this figure shows the effect of temperature on both 
the measured and the pairwise average Xxz- 

V. CONCLUSIONS 

We have developed a spin-pair model for understanding the behavior of dilute 
LiHo :r Yi_ a ;F4. A weighted ensemble average of all spin pairs reproduces the high-transverse- 
field experimental susceptibility, but not the low-field antiferromagnetic character of the 
data. Nonetheless, the rise in the longitudinal susceptibility at a transverse field of around 
1 T, which looks like a signature of a spin gap, does correspond to the calculated suscep- 
tibility for certain antiferromagnetic pairs. This suggests that a full understanding of the 
system requires treatment of larger clusters, an extension which should be numerically fea- 

15 



o 1 ^ 

0.0 0.5 



0.5 1 1.5 2 0.5 1 1.5 2 

Transverse Field (T) 



FIG. 9: (Color online) Measured and computed Xzz (in units of emu/mol Ho), (left) Computed 
(heavy, red curves) and measured (blue, light) susceptibility. Solid, dotted and dashed curves are 
T=70 mK, 110 mK and 150 mK respectively. (Right) The effect of tilting the transverse field 
(dotted curve), and of omitting the hyperfine interaction (dashed), at T=70 mK. Unbroken red 
(heavy) and blue (light) curves, again show the computed (with hyperfine, no tilting) and measured 
susceptibilities respectively. 

sible because of the observation that the primary effect of the hyperfine splitting in the dc 
susceptibility is to renormalize the transverse field. This will allow extension of the model 
to larger clusters of spins using the simplified 2-state description for individual spins rather 
than a full 16-state description. Ultimately, to reach the thermodynamic limit, it would still 
be necessary to generalize a scaling approach, such as the real-space renormalization group 
of Ref. 9 , to include finite transverse fields. 

Such an extension would sample somewhat different regions of configuration space, since 
Fig. [7] shows that the antiferromagnetic region extends considerably farther in distance than 
does the ferromagnetic region. This space is not sampled significantly in the pairwise model, 
owing to the rapid fall-off of the weighting function Wi with distance, but larger clusters can 
sample this interaction region far more extensively. 
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FIG. 10: (Color online) The apparent transverse susceptibility resulting from a transverse field 
tilted by 0.6°. (Left) Measured susceptibility (blue, thin, solid curve) is contrasted with the sus- 
ceptibility of a single ion in a tilted field (dashed curves), and the pairwise average susceptibility 
(dotted). The heavy, red curves include hyperfine effects. The thin, red curves do not. (Right) 
The effect of temperature. Measurements are shown as light, blue curves; calculations as heavy, 
read curves. Temperatures are 70 mK (solid curve), 110 mK (dotted) and 150 mK (dashed). 
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